/*
Network Design Project
"Minimum Bound Degree Spanning Tree"
Version: 0.4
Authors: Alessandro Baroni, Paolo Bozzelli, Giuseppe Cascavilla
Domain: univaq.it
*/

#include <stdio.h>
#include <stdlib.h>
#include <string.h>
#include <glpk.h>
#include "lib/ourlib.h"

glp_graph *G;
int n,e,last_sol_allocated,starting_solution_allocated,arrays_len;
double *starting_solution;
double *last_sol;

int solve(int);
int	diff(double*,double*);

int main(void) {
	char buffer_str[128];
	const char* GRAPH_FILE = "generated_graph.txt";

	int i,is_changed;
	
	G = custom_create_graph();
	i = custom_read_graph(G,GRAPH_FILE);
	e = G->na;
	n = G->nv;
	if(i != 0){
		//TODO error handling
		return -1;
	}
	
	last_sol_allocated = 0;
	starting_solution_allocated = 0;

	arrays_len = 0;

	for(i = 1; i <= n; i++){
		//printf("risolvo al nodo %d\n",i);
		//getchar();
		solve(i);
	}
	
	
	if(arrays_len != 0){
		starting_solution = malloc(sizeof(double)*arrays_len);
		for(i=0; i<arrays_len; i++){
			starting_solution[i] = last_sol[i];
		}	
	}
	
	is_changed = 1;
	while(arrays_len != 0 && is_changed){
		printf("Eseguo di nuovo\n");
		fflush(stdout);
		system("./separazione.out");
		//getchar();
		
		free(last_sol);
		last_sol_allocated = 0;
		
		for(i = 1; i <= n; i++){
			//printf("risolvo al nodo %d\n",i);
			//getchar();
			solve(i);
		}
		
		if(last_sol_allocated == 0){
			//se non ho trovato nessun subtour
			is_changed = 0;
		}else{
			//se ho trovato almeno un subtour
			is_changed = diff(last_sol,starting_solution);		
			for(i=0; i<arrays_len; i++){
				starting_solution[i] = last_sol[i];
			}
		}

	}
	free(starting_solution);
	
	printf("La soluzione non cambia\n");
	
	return 0;
}

int diff(double*x ,double* y){
	int i;
	//printf("arrays_len = %d\n",arrays_len);
	for(i=0; i<arrays_len; i++){
		printf("%d %lf %lf\n",i,x[i], y[i]);
		fflush(stdout);
		if(x[i] != y[i]){
			return 1;
		}
	}
	return 0;
}

int solve(int num){
	glp_prob *lp;
	glp_prob *separazione;
	const char* SEPARAZIONE_TEMPLATE = "separazione%d.lp";
	const char* CONSTRAINT_TEMPLATE = "subtour elimination constraint for z_%d - %d";
	const char* LP_FILE = "mbdst-formulation.lp";
	const char* SOLUZIONE_FILE = "soluzione";
	char varname[15],separazione_file[128],constraint_name[128];
	double z,*x,*coeff_array,rhs;
	int i,j,rownum,subtour_cardinality,*indexes_array;
	
	strcpy(varname,"");
	sprintf(separazione_file,SEPARAZIONE_TEMPLATE,num);	

	
	lp = glp_create_prob();
	i = glp_read_lp(lp,NULL,LP_FILE);
	if(i != 0){
		//TODO error handling
			return -1;
	} 
	glp_create_index(lp);
	
	separazione = glp_create_prob();
	i = glp_read_lp(separazione,NULL,separazione_file);
	if(i != 0){
		//TODO error handling
		return -1;
	}
	
	glp_simplex(separazione, NULL);
	
	z = glp_get_obj_val(separazione);
	

	x = malloc(sizeof(double)*glp_get_num_cols(separazione));
		
	for (i = 0; i < glp_get_num_cols(separazione); i++) {
		x[i] = glp_get_col_prim(separazione, i+1);
		
	}
	
	
	subtour_cardinality = 0;
	
	for (i = 0; i < glp_get_num_cols(separazione); i++){
		strcpy(varname,glp_get_col_name(separazione, i+1));
		if(x[i] == 1.0 && varname[0] == 'w'){
			subtour_cardinality += 1;
		}
	}
	
	if(subtour_cardinality >= 1){
		rownum = glp_get_num_rows(lp) + 1;
		glp_add_rows(lp,1);
		sprintf(constraint_name,CONSTRAINT_TEMPLATE,num,rownum);	
		glp_set_row_name(lp, rownum, constraint_name);
		glp_set_row_bnds(lp, rownum, GLP_UP, 0.0, subtour_cardinality - 1);	// <= subtour_cardinality - 1
		indexes_array = malloc(sizeof(int)*(subtour_cardinality + 1));
		coeff_array = malloc(sizeof(double)*(subtour_cardinality + 1));
	
		
		// Stampa su stdout del risultato
		//printf("SEPARAZIONE\n");
		//printf("\nObjective :\tz = %g;\n", z);
		//printf("Vector:\n");
	
		j = 1;
		for (i = 0; i < glp_get_num_cols(separazione); i++) {
			//printf("\tx%d = %g;\n", i+1, x[i]);
			//printf("\t\t%s = %g\n", glp_get_col_name(separazione, i+1), x[i]);
			strcpy(varname,glp_get_col_name(separazione, i+1));
		
			if(x[i] == 1.0 && varname[0] == 'w'){
				// x[i] == 1 nel problema di separazione
				//la variabile associata alla colonna i+1 è nel s.e.c. da aggiungere a lp
				varname[0] = 'x';
				//printf("add: %s\n",varname);
				indexes_array[j] = glp_find_col(lp,varname);
				coeff_array[j] = 1.0;
				j++;			
			}
		}
	
		glp_set_mat_row(lp,rownum,subtour_cardinality,indexes_array,coeff_array);
		
		glp_simplex(lp, NULL);
		
		// qui aggiornavo la starting_solution
		
		if(last_sol_allocated == 0){
			arrays_len = glp_get_num_cols(lp);
			last_sol = malloc(sizeof(double)*arrays_len);
			last_sol_allocated = 1;
		}
		
		z = glp_get_obj_val(lp);
		// Stampa su stdout del risultato
		printf("PROBLEMA\n");
		printf("\nObjective:\tz = %g;\n", z);
		//printf("Vector:\n");
		for (i = 0; i < glp_get_num_cols(lp); i++){
			//printf("\tx%d = %g;\n", i+1, x[i]);
			//printf("\t\t%s = %g\n", glp_get_col_name(lp, i+1), glp_get_col_prim(lp, i+1));
			last_sol[i] = glp_get_col_prim(lp, i+1);
			strcpy(varname,glp_get_col_name(lp, i+1));
		}
		printf("Aggiorno last_sol\n");
		
		glp_write_lp(lp,NULL,LP_FILE);
		custom_write_sol(G,lp,SOLUZIONE_FILE);
		
		free(indexes_array);
		free(coeff_array);
	}
	
	glp_delete_prob(lp);
	glp_delete_prob(separazione);
	
	return 0;
}
